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I, INTRODUCTION AND SUMMARY 



A. BACKGROUND 

The Kslmsn filter is a recursive Beyesien lesst-souares 
estimator of an n-dimensional system state vector based on an 
m-dimensional measurement vector. The filter may operate in 
a J-dimensional coordinate system where J:i;;m !» Jiiin. The basic 
assumption is that each dimension of the coordinate system 
varies according to a kth order Gauss-Markov process. The 
Kalman Filter was developed in the early 1960's by Kalman and 
Bucy Crefs. 1 and 23. 

/ 

The Kalman filter may be used to obtain an optimal 
estimate of the present statey a prediction of future states; 
and/or smoothed estimates of past states. The current state 
estimate is Generally used to determine an optimal control 
input. Future state estimates are used to determine optimum 
present policy. Smoothed past state estimates are used for 
data analysis and model buildind. Thus the potential areas of 
application span the field of time series analysis. 

Applications of the Kalman filter are numerous and the 
theory is beind continually developed and extended. An 
overview of the development of linear filterind theory and an 
extensive bibliodraphy may be found in Kailath C ref. 33. A 
reasonably clear presentation of theory and applications is 
contained in Gelb Cref.43. 
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Perhsps the widest and most successful application of 
Kalman filtering has been to vehicle tracking and control* 
Clark C ref *5.1 has written a particularly lucid description of 
the desidn of a filter for an ant i~ai rcraf t dun fire control 
system which is noteworthy for its clarity of presentation of 
the underlyind theory. It is evident Crefs. A and 5U that 
the desidn process is heuristic j and reouires extensive 
testind and analysis of candidate filter conf idurations? even 
when the process is well-understood and is based on a mature 
technolody ♦ 

The Kalman filter has also been appliedK with varyind 
dedrees of success y in economic models x inventory models i- and 
even weather models. Considerable difficulty is encountered 
in model buildindr because the filter desidn reeuires dood 
estimates of the variance and covariance of noise sources f as 
well as an accurate state transition model. A prior estimate 
of system state and covariance is also reouiredr which is 
somewhat less critical because errors in the prior estimate 
decrease with time. These parameters are often difficult to 
determine in hidhly random processes of Questionable 
stat ionar i ty . 

The Kalman filter is derived and desidned almost entirely 
within the time domain^ althoudh Clark Cref.53 does refer to 
the concept of filter bandwidth. The Kalman filter is 
essentially a low-pass filter with a very wide transition 
bandf and hidher-order filters have some ampl if ication at the 
mid or low-mid freauency rande. In deneralf the stop band 
does not completely attenuate hidh freauencies. This allows 
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the filter to otteriuate hish-f reouency noise somewhat while 



still retaining some response to sudden changes of state. 



B. PURPOSE 



The purpose of this thesis is to acouaint the reader with 
the Kalman filter? to show how the choice of various filter 
parameters affect its performance? and to provide design 
insisht through analysis in the freouency domain. The 
approach is tutorial? and the reader is referred to some of 
the interesting examples which may be found in the 
literature . 

C. METHOD 



The freouency response of several simple filter designs 
were investigated usins the Fast Fourier Transform prodram in 
the APL Library 2. The computer results were Justified 
analytically for the simplest desidn? a scalar sindle-state 
filter. Derivations are presented in appendix A. 



D. LEUEL OF PRESENTATION 



Full understandi nd of the theory reouires 
stochastic processes that evolve over time? 
understandind of didital sidnal theory in 
domain. The Fourier transform is a basic 



a knowledde of 
as well as an 
the freouency 
tool. A full 



exposition of the under lyind theory is clearly beyond the 



scope of this presentation. The reader is directed to Larson 
and Shubert Cref*6J for the theory of stochastic processes 
and to HamminS Cref.7H for the theory of digital filtering. 
As previously mentioned^ Gelb Cref.43 and Clark Creft5D are 
Sood references for the Kalman filter. Bloomfield Cref.83 

and Brillinder Cref.91 are also applicable references. Brown 
C ref. ion and Box and Jenkins Cref.llU contain related 
material . 

There are few readers who are entirely conversant with 
both the freouency domain and time domain approach to time 
series analysis. Nevertheless s' a duality exists between the 
two> and a summary of the theory is presented. 

Illustrative examples will often be based on trackind 
models» because this is presently the widest area of 
application of Kalman filters? and because most readers will 
find the concepts of position? velocity? and acceleration 
easy to understand. The concepts are easily extendable to 
other areas. For example? the economist may wish to replace 
"velocity" with "trend". 

E. SUMMARY OF RESULTS 

The steady-state dain? bandwidth? and sensitivity of the 
linear discrete Kalman filter are shown to be completely 
determined by the choice of the process and measurement noise 
covariances. Filter performance on stationary or 
nearly-stationary data can be predicted by comparind the 
freouency response of the proposed filter with a spectral 
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analysis of the data* The wide transition band of the 
amplitude response of the scalar Kalman filter can be 
sharp'ened by multiple passes of the data through a 
hidher-dain filter* This can be accomplished simply and 
recursively. The superiority of symmetric smoothind filters 

over non-symmet r ic filters was demonstrated. When used as a 
smoother <by usind both forward and backward passes) the 
Kalman filter was as effective as a non-recursive Gaussian 
filter. Hidhet — order filters were shown to have hidher 
bandwidth and ampl if ication as the order of the filter was 
increased. A freouency domain approach to filter desidn may 
provide additional insidht and enable the desidner to achieve 
better filter performancey particularly when the system state 
transition model and noise covariance models are not well 
understood , 
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II 



THEORY 



A. STOCHASTIC PROCESSES AND STATIDNARITY 

A continuous stochastic process X(t) is 3 Osussian 
process if the probability densities of all orders are 
multivariate Gaussian densities. It is a kth order 
Gauss-Markov process if the state at time t depends only on k 
earlier states. If we should expand the state space to k 
states^ which include all derivatives up to the <k"l)th7 the 
future system state vector will depend only on the present 
state. For examplof if the acceleration of a vehicle is a 
first-order Gauss-Markov p-rocessi’ then the position of the 
vehicle is a third-order Gauss-Markov process. Howeveri- if 
our state space includes acceleration and velocity as well as 
position? the future state of the system is independent of 
all but the present state. If the random acceleration has 
zero mean? and variance one over one time increment? the 
acceleration is a standard Wiener process W(t). The 
derivative of the Wiener process? written dW(t)? has zero 
mean? unit variance? and is called white Gaussian noise? 
which may be thought of as a *zero~th order* Gauss-Markov 
process Cref.43. 

The standard Wiener process is not stationary? because 
the variance drows linearly with time. That is? the estimate 
of a future state based on the present state has variance 
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that is 3 linear function of time. However > the standard 
Wiener process hes stationaryi- independent increments. That 
is? the variance at time <til) diven the state at time (t) is 
constant and independent of t. 

A stochastic process X(t) is wide-sense stationary if and 
only if it has a constant mean function? and a correlation 
function such that C ref. 6.1 

R, (t, + is? t^ f s) = R (t, ? t, ) = R (t^-t, ) 



that is? the correlation function of the process is 
independent of an arbitrary time shift s. A Gaussian process 
is strictly stationary if and only if it is wide-sense 
stationary Cref.61. 

The Gauss-Markov assumption makes possible the 
development of theory and appl icat ions ? because? in deneral? 
any linear operation performed on a Gaussian process results 
in another Gaussian process? and the Markov property allows 
consideration of only the present state? disredardind all 
previous states. 



B. THE PHILOSOPHICAL CONCEPT OF STATIONARITY 

A f reauency-domain analysis of a stochastic process is 
only meanindful if the process is stationary. If the process 
were chandind over time? the spectrum would chande over time. 
Since the spectrum can only be analyzed by means of data 
taken over time? such analysis of a non-stationary process 
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would be iiieenindless ♦ However v if the process is 
"euesi-stetionary * » that is;* it exhibits stetioners 
statistics for s while)* then underdoes 3 chsnde* then settles 
down to stationsri tvj adain» the freeuencs approach is still 
useful > althoudh inaccurate over the transition period. As 
an example? consider an airplane subject to random 
accelerations due to air turbulence. An appropriate model 
midht be a third-order Gauss-harkov process as lond as the 
airplane maintains a straight path or turns at a constant 
accelerat ion . However? the pilot's inputs to initiate or 

terminate a maneuver would result in brief periods of 
non-stationarity ? and the model would perform inadeouately 
durind and immediately after the transition period. 

It may be ardued that every practical process can be 
considered stationary over infinite time. If the process is 
random? it represents an ensemble of possible paths? of which 
any realization in terms of real-world data is only one 
possible Path? and may or may not be closely representati ve 
of the ensemble. When dealind with reality? we are often 
forced to assume stationarity in order to make analysis 
possible? and often we obtain dood results even thoudh we can 
never know whether or not the assumption of stationarity is 
really valid. 

C. DUALITY OF THE TIME AND FREQUENCY DOMAINS 

1 . Fourier Series 

A very wide class of mathematical functions may be 
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rer-- resented by the Fourier series II ref. 123 3S follows^ 
d(t) ••= 3^ + ^_(^^cos nt f b sin nt> 

h»< 

where 0 :i t ii! 2 7T. 

E;;istence <3nd convergence of this series reeresentstion 
reouire only that <2<t) be everywhere sindle-valued!- and 
possess a finite number of maKiiiiS!' minima!' and finite 
discontinuities. The function d(t) need not be 
differentiable. Any function meetind the above criteria can 
be thoudht of as a constant mean function a^j plus an 
infinite series of sines and cosines of intedral freouencies 
and various amplitudes. Of courser the independent variable 
t must be shifted and scaled to the interval C 0r2rT3. Note 
that the lowest freouency present? aside from the 
2ero~f reouency mean? is one cycle for the span of d<t). 
Amond the functions mectind the criteria are a souare pulse? 
an impulse? and any manif estation of a random walk. In 

practice? the Fourier analysis of a function d(t) reouires 
the truncation of the infinite Fourier series. This results 
in a smooth least-souares approximation to the function d(t). 
There are ripples in the approximat ion if the function d(t) 
is not differentiable or if the truncation is too severe. 
This is known as the Gibbs phenomenon? and is illustrated in 
fidure 1? which was taken from Hammind II ref. 73. By takind a 
sufficient number of terms in the Fourier expansion? we can 
improve the closeness of the approximat ion . 
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2 . Basic Concept 

The basic concept of the duality between the time and 
freouency domains is so simple that it often Sets lost in a 
forest of Fourier transforms* The time period is the 
reciprocal of freauency* The basic relationship is 

V/27T = f = 1/T 

where v is the freouency in radians/unit time? f is the 
freouency in cycles/unit timer and T is the time period for 
one cycle* Stated simplyr freouency is the inverse of the 
time period* 

3 * Discrete Data and the SamplinS Theorem 

The disital computer allows the efficient analysis of 
continuous phenomena by means of discrete approKimations * We 
saw earlier that the lowest freouency contained in a Fourier 
expansion of a function S<t) was the reciprocal of the time 
span covered by the function* Similarlyr the famous Samplind 
Theorem IZrefs* 6 and 73 states that if a function S(t) in 
continuous time is sampled at constantr discrete time 
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int©rv<3ls (that isy at a rate of ,1/tvt)!' then the highest 



observable freauency is 0»5 cycles per nrieasurement interval 
Lt. This means that at least two observations are reouired in 
each cycle in order to observe that particular freouency. 
The freouency 0»5/tty usually written simply 0,5y is 

referred to as the Nyauist freeuency* The result is the 
aliasing phenomenony which is familiar to most moviegoers* 
Durind the chasey the stadecoach wheels appear to stop or 
rotate slowly backwards when the rate of rotation of the 
wheel spokes (spokes/sec) exceeds 1/2 the camera rate 
(frames/sec)* When hidher freoiuencies exist in the function 
d(t) sampled at a rate t>ty they are folded back and appear in 
the freeuency spectrum of the sampled data as freeuencies 
less than the Nyouist freouency. The samplind theorem shows 
that a spectral analysis of discrete data is only meanindful 
over the Nyauist interval C -0.5/Iity fO*5/LtII* 

A ♦ The Discrete Fourier Transform 

Any function d(t) for which a converdent f’ourier 
series exists may be represented in the freauency domain in 
terms of real and imadinary parts y or in terms of amplitude 
and Phase andley as a function of freauency* It should be 
noted that the function d(t) may also be complex-valuedy but 
we will deal with only real-valued functions* In the 
continuous domainy the formulas 

«o 

d(t) = irf-J" exp(ivt) dv 
and ' 
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G<v) 



dt 






= fd(t) 



exp ( "i V b ) 



represent 3 Faurier trsnsform pair» The frequency response 
G(v) completely determines the time function S(t) and 
conversely . 

If the function d(t) is sampled at intervals t = 
0 » 1 » 2 i' . ♦ . »n the time-to-f reouency transf ormation becomes 

/I 

G<v) =l^d(t) exp < - i vt ) t - 0 y 1 ? ♦ . ♦ » n 

t*o 

which is defined only on the Nyauist interval L~TT r TT 3 , here 
defined in radians* The formula may be written in a more 
familiar form by usind the Euler relation 

exp(-ivt) = cos vt - i sin vt 



as 



n 

G(v) ~ d ( t ) (cos vt ” i sin vt) 

t=o 



which is c 


ontinuous in v on the 


interval 


C-TTfiri, The 


Fourier 


transform 


is a hit difficult to 


handle 


analytical ly 


for all 


but the 


simplest functions!' 


but the 


discrete 


Fourier 


transform 


is denerally easy to 


compute 


by use of 


a Fast 



Fourier Transform (FFT) prodram available in most computer 
libraries* The output will denerally be a very close discrete 
approximat ion to G(v)r if the span of d(t) is larde enoudh* 
The inverse t ransf ormation can also be made* 
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Since G(v) is complex valued whenever the function d(t) 
is not symmetricr it is often useful to represent it in terms 
of amplitude and phase. The amplitude is 



|G<v)| = J 



G(v)G(-v) 



/ 



CReCv)/ 



f 






where Re(v) and Im(v) are the real and imaginary parts of 
G(v)« The phase andle is 



0(v) = arc tan Im(v) / Re(v) 



D, THE DOOB-MEYER-FISK DECOMPOSITION 



In most practical applicationsy 
sample-continuous stochastic process X(t) 



a f inite-variance 



c a n be w r i 1 1 e n 



X(T) 



X(0) 



T 

t f Jl 



r 

B ( t ) dW ( t ) 



d c 

where X(0) is the initial value of that processy A(t)dt is 
predictabley smooth behavior determined by a set of 
deterministic differential eouations describing the system y 
and B(t)dW(t) is noisey where dW(t) is white Gaussian noisey 
and B(t) is a smooth t ransf ormat ion that is sometimes thought 
of as "coloring* the noise. Such a representation is called 
the Doob-Meyer-Fisk decomposi tion Cref,6Iy which may be 
thought of as separating the process into a sidnal and noise. 
Several important points must be made with redard to this 
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GQij3tion» It is not intended here that the eicF'ression he 
evaluated analytically* The integral B<t)dW(t) is an I to 
integral .* which is not even a stochastic version of a 
StieltJes intej^ral Cref*6J* Also» slthoudh the processes 
A(t) and B(t) are smooth functions that may be considered 
deterministic representations of system behavior^ they are 
not necessarily known to the observer? even when an adeouate 
technological representation exists* 

Consider adain our piloted aircraft beind tracked by a 
radar* The process A(t) represents the dynamics of the 
airframe? as affected by the control inputs of the pilot? 
which are unknown to the radar observer* The process B(t) 
consists of several parts* One is the measurement process? 
which may or may not be known to the radar observer* For 
example? a nutatind radar antenna midht impose some periodic 
error in the measurement? which would be manifested in the 
process EUt)* Overlaid on this midht be a white Gaussian 
noise measurement error* Air turbulence could also be 
represented as white Gaussian noise? which? however? could 
only be manifested throudh deterministic airplane dynamics* 
There are those who would ardue that the pilot should also be 
modelled as a random variable* In any event? the process 
E<(t) midht be further decomposed into several processes? here 
at least airframe response to air turbulence and periodic 
radar antenna dv^namics* 

The vital observation is that if the freauency content of 
the processes A(t) and E<<t) are known to be different? they 
can he Partially separated by a spectral analysis of the 
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Mi 





•f reauencies in «1.1. control axes 



These can he estimated 



closelyy even for enear^j airFlanes? and are scnerally similar 
amond similar types of airplanes.* althoudh they vary with 
airspeed. .It is physically impossible for the airplane to 
respond faster than its hidhest natural dynamic freouencies. 
Any freciuency content hidher than this must he noise. If the 
radar system dynamics are of a hidher freouency than this? 
they can also he separated. The pilot will take advantade of 
the full response rate of the airplane only very rarely. 
Therefore!- low freeuency components are most likely due to 

Pilot maneuvers. Of courses- since white Gaussian noise has 
a flat freeuency spectrum as a result of aliasind C ref .73* it 
is impossible to separate all of the noise from the sidnal. 
However)- it is often possible to remove auite a bit of it. 

E. DIGITAL FILTERS 

A didital filter is a linear t ransformation applied 
iteratively to a set of data points. The purpose here is to 
separate noise from the sidnal. The simplest didital filter 
is the simple averade)- which estimates the mean value from 
the data I- and smooths out all fluctuations. The most deneral 
form of the didital filter was stated by Hammind Cref.73 as 



where the estimate x<t) at some point t is a linear 




22 



combi nation of the data nointa z(t-k)i> and r-erhans of the 
previous estimates x(t~k)* The coefficients a(k) and b(k) 
are weidhtind coefficients and may? of course i> be zero. As a 
result of the samplind theoremy the filtering process is 
meaningless unless the measurements z(t-k) are made at 
eoually spaced intervals alond the t axisr where t is 
usually? but not necessarily? time. 

1 . Some Classif icat ions of Didital Filters 

Didital filters may be classified as symmetric or 
non-svimmetric? and as recursive or non-recursive. A 
symmetric non-recursive filter is one in which all b(k) eeual 
zero and all a(k) = a(-k)? such as the filter 

x(t) = 0.2 z(t-l) + 0.6 z(t) f 0.2 z(tfl). 

An example of a recursive filter is 

x(t) == a z(t) ■{■ b x(t-l) 0<3<1? b = 1-a 

which is not symmetric. This particular filter may be 
expressed as 

x(t) = a z<t) f bC a z(t-l) + bL' a z(t-2) + ....3 3 
which reduces to 

x(t) •- a z(t) + ab z(t-l) + ab^z(t-2) +...+ ab’'z(t-n) i . . . 
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The recursive filter extends to the infinite P3st» slthoudh 
the coefficients 3b will opprosch zero? if |h|<l» In this 
cose? 3 recursive filter csn he closely spproximsted by 3 
non-recursive filter* A primary advsntsde of the recursive 
filter is that old data need not he stored. New estimates 
may be computed simply and rapidly as time evolves. This is 
an important advantade for real-time 3pp1 ications . 

2 . Applications of Didital Filters 

Didital filters are used to separate a sidnal from 
noise? to separate various freeuency components of a sidnal? 
and/or to perform such mathematical functions as intedration 
and differentiation. A review of Simpson's rule and the 
Trapezoidal rule should convince the reader that these 
numerical intedTation technieues are? in fact? recursive 
didital filters. Sometimes a filter has two purposes. For 

I 

example? it midht be desirable? in estimatind velocity from 
successive observations of position? to simultaneously 
differentiate and remove hidh freouency noise. Uhen a filter 
is used to stop part of the freouency spectrum? it is 
referred to as a "low-pass"? "hidh-pass"? "band-pass"? or 
"band-stop" filter? dependind on its function. 

3 . Analysis of Didital Filters 

In the time domain? a didital filter is described 
completely by its impulse response function? which is nothind 
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more than the response of the filter to deta consisting of a 



strings of :ceros)' a sindle one)- followed by zeros. The output 
of the filter is then simply the weidhtind coefficients a(k). 
If the filter is recursive y we midht not be able to deduce 
the recursive form from the coefficients a<k)» but that will 
not concern us here. The Fourier transform of the 
i m p u 1 s e - r e s p o n s e function 

N 

H < V ) = a ( k ) exp ( - i vk ) 

k=-H 

will completely specify the freouency response of the filter. 
If the filter is symmetricy there will be no imadinary party 
and hence no phase shift. If the filter is recursive y it 
cannot practically he symmetric y and the summation will 
denerally run from zero to infinity. That isy the impulse 
response will e>rtend infinitely far into the future y which 
means that the filter remembers all of the past. 

The duality of the time and freouency domains allows us 
to specify a desired freouency response and to desidn an 
appropriate filter by calculatind filter weidhtsy or to 
analyze an existind filter by calculatind the freouency 
response from the filter weidhtind coefficients. 



F. DATA ANALYSIS AND EXPERIMENT DESIGN 



No didital filter should be applied to data 
without a clear idea of the effect of the filter 
data. Slutskvj and Yule first noted that some 



analysis 



upon the 
smoothind 
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formulas induced ?'eriodic functions 



in the smoothed estimate 



that were more the effect of smoothing than of the original 
data C reft 73 ♦ A spectral analysis of representative raw data 
can be helpful in deciding on an appropriate filtering 
technique » However? such data as economic time series or 
weather data are typically very noisy? are based on a 
relatively short run of data? and cannot be described by an 
adeouate technological model* The analyst must be aware of 
these problems. Sometimes there are no dood solutions? but a 
spectral representation may produce freouencies that can be 
explained on rational drounds. 

Another potential pitfall is a result of the samplind 
theorem. Consider the timely example of an air pollution 
model. It would he reaonable to suspect that air pollution 
would follow at least a daily cycle? or perhaps an eiriht hour 
cycle if mornind and evenind rush hours were considered. 
Daily samples of air pollution could not hope to uncover 
cycles of a shorter period than every two days. Samples 
every four hours would be mardinally adeouate. Hourly 
samples would he necessary for a dood analysis. Additionally? 
recall the reouirement for eoual ly-spaced samplind intervals. 
For various reasons? the analyst may have no control over 
data collection. However? he must always understand what has 
been done? or could have been done? to the data? as well as 
what he is doind to it? in order to avoid erroneous 
conclusions . 
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THE LINEAR DISCRETE KALMAN FILTER 



A. DESCRIPTION 



The linear discrete Kalman filter is a recursive Bayesian 
least-sQuares estimator of the state vector of a linear 
system based on a vector of noisy measurements made at 
discrete time intervals. The process to be estimated is 
assumed to be an n -state Gauss-Markov process of order k? 
subject to process noise U with zero mean and covariance 
matri;-; Q. The process is observed by an m-dimensional 
measurements subject to measurement noise V (not to be 
confused with freeuency (v>) with zero mean and covariance 
matrix R» The filter reouires a prior Bayesian estimate of 
system state and covariance. The recursive estimate of system 
state at time t is obtained by the formula > 



X(t!t)= X(t!t--1) + K(t)CZ(t)-HX(t!t--l)l 



where 



X(t! t) 
Z(t) 

X( 1 1 t“l ) 
K(t) 

H 



state estimate based on current measurement 
current measurement 

state estimate prior to current measurement 
Kalman dain matrix (to be discussed later) 
observation matrixy which is constant 



The derivation of the Kalman filter eeuations may he 
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found in Gelb C ref* -411. A suirimsry of the filter eaustions is 
presented in fidure 2? which should be consulted in order to 
follow the subseeuent discussion* 

In denersly the st<3te model represents 3 dynamic ssstemy 
that isy one which chandes with time* The extrapolation of 
the state estimate to the time of the next observation is 
obtained by the formula J 

X<t-M It) = § X(t! t) 

where ^ is the state transition matrix* The observation 
process occurs accordind to the conceptual relation 

Z(t) == HX(t) + K> 

where X(t) is the true system state observed throudh the 
observation matrix Hy and V represents measurement noisey 
which is assumed to be a Gaussian random variable with zero 
mean and covariance matrix R* Note that the process 
represented by this formula is assumed to occur in the real 
world* The computation does not occur in the filter* Rathery 
the measurement Z(t) is an input to the filter* 

In the linear Kalman filtery the da in K<t) does not 
depend in any way on the data* It depends only on the model y 
and is therefore extremely sensitive to assumptions* Gain is 
calculated accordind to the formula 

T T 

K(t) = P<t)H CHP(t)H + R3 
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PiG-uae- 




where P(t) is the covariance in the system state estimate 
prior to the current measurement and R is the covariance of 
the measurement error* The covariance is updated according 
to the formula 



r(t) = CI-K(t)H3 P(t) 



where S(t) is the state covariance diven the current 
measurements and I is the identity matrix* The covariance is 
extrapolated to the time immediately prior to the next 
observation by the formula 

P<tfl ) = ^ 2.<t > + Q 

where ^ is the state transition matrix and Q is the 
covariance of the process noise* Combining the above two 
eauations shows that the covariance of the state estimate at 
the time of the current measurement depends on the previous 
covariance accordind to the formula 

Z(t) = i;i-K<t)H.l C(J r<t-l) f G}.-] 



Filter performance is very dependent on adeauate 
modellinds particularly on the state transition model ^ and 
the choice of noise covariances R and Q* To a lesser extents 
performance also depends on the initial estimates of system 
state X(1!0> and covariance P(IIO)* Howevers the latter 
parameters are less important because their effects decrease 
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with time 



If the matrices Rs' Qy $y 



and H are constant 



in 

timey the i^ain K(t) and covariance matrices S<t) and P(t) 
eventually reach a steady statey and are completely 
determined by Ry Qy ^y and H* 

For a 5iiven linear filtery it will be shown that filter 
slain y covariance y and freauency response will be completely 
determined by the choice of R and Q» 



B. THE SCALAR KALHAN FILTER 



The multi-state Kalman filter 
device* Howevery it is difficul 
manipulate in closed form because 
of singular matrices. An analysis 
filter can be used to illustrate 
filtery and to aid in developing 
In the discussion that followsy 
matrices are scalarsy andy in part 
Matrix notation is preserved for c 
found in appendix A* 



is a powerful computational 
t and often impossible to 
of the freouent occurence 
of a sinsle-state (scalar) 
the mechanics of the Kalman 
an intuitive understandinS ♦ 
it is assumed that all 
iculary Q and H eoual one. 
larity. Derivations may be 



1 . 



dain can 



Transient and Steady-State Gain 
It can be shown (appendix A) that 
be expressed recursively as 



the scalar Kalman 



K(t) = K(t-l) + Q/R 

K(t-l) + Q/R f 1 

Uhen the filter reaches steady-state y the Sain is constant 
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The inverse relstionship is 



£ 

R 



K 

1-K 



Thusi» the variance ratio Q/Ry which is the ratio of process 
noise variance to nrieasu rement noise variancey completely 
determines the steady-state Sain. The steady-state filter is 
completely described by the formula 



X(t) = K Z(t) f <1-K) X(t-l) 



2 » Frequency Response 

Lettind K = a and <1-K) •- by the impulse-response 

function G(t) may be written 



•t 

G(t) = ab y t=0 y 1 y 2 y ♦ * ♦ ♦ ♦ 
The Fourier Transform is 

H(v) = JciCt) e:;p (-ivt)dt 

H<v) ••= 3 ^Cb exp<-iv)3^ 

tso 



H( v) 



a / i:;l 



b e>;p ( - i V ) 3 



Since the filter is not synrimetrici' the freauency response 
H(v) has both resl and iinodinary parts* The amplitude may be 
written 

A = |H(v) I 
which reduces to 



The Phase andle may be written 

0(v) - arctan / -b sin v | 

\1 -b cos V/ 

The andle for maKimum phase shift is 



= JH ( V ) H ( -V ) = a / J~i + b^ ~2b cos v 



A == / Q/R 

J Q/R f 2<1 - 



cos v) 



v(ma>{ 0) = arccos b = arccos (1-K) 



V (max 0) = arccos (1 f Q/2R 



yQ*V4^ 



f Q/R ) 



Thereforer the variance ratio Q/R also completely specifies 
the steady-state freeuency response of the filter* Amplitude 
and phase relationship for several values of dain are plotted 
in figures 3 and 4* 

It is evident that hidh Q/R (hidh dain) reduces the phase 
lad of the filter but allows more of the hidh~f reouency 
components to pass* Conversely^ low Q/R (low dain) 
attenuates more of the hidh-f reouency components >« at the 
expense of an increased phase lad* Note that even at very 
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low .^ain (low Q/R) not all of the hi<2h-f reouency component is 



attenuated? and the phase lad is Quite severe. The slope of 
the amplitude chande is Quite shallow? imply ind that 
attenuation increases dradually as freauency increases. This 
is 3 conseouence of the assumptions and performance will not 

he adeauate if the data does not represent a Gauss-ha rkov 
process? but in fact represents some phenomena chandind with 
t ime . 

C. IhRROVING THE FILTER 

The transition band of the filter can be sharpened? and 
more of the hidh- freauency components eliminated? by runnind 
the data throudh two filters in series. The basic scalar 
filter was 



x(t) = a z(t) + b K(t~l) 

where a == K and b = (1-K). Runnind the data throusth the 
filter adain? we obtain a new estimate y(t) ? where 

y(t) =3 x(t) lb y(t~l) 

It should be evident that we can accomplish this all in one 
step as 

y(t) - a^z(t) f 2b y(t-l) - b y(t-2) 
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We need only to save one additional previous estimate y(t~2) 
as well as y(t-l). The impulse response function is 



d<t) = (tfl)a* » t = Oi-lr... 



We have performed a convolution in the time domain? which 
corresponds to a mul tipi ication in the freouency domain* 
This may not he e;;actly what we want* Let us suppose that we 
want the weidhtind coefficient for the present data point 
z(t) to be 0*27 in both cases* This reeuires a = 0*27 for 
the basic filter and a = ^0 * 27 - 0*52 for the double filter* 
The impulse response function for both filters is presented 
in fidure 5* Note that the double filter fordets the past 
more readily* The amplitude and phase shift for both filters 
is presented in fidures 6 and 7* The dain for the scalar 
filter was 0*27» correspondind to a variance ratio ( Q/R ) of 
0*1* Also note that the double filter has somewhat better 
hidh-f reouency attenuation? somewhat less attenuation at low 
freauencies? and a slidhtly sharper (steeper) transition* We 
would therefore expect it to be a hit better at separatind a 
low~f reouency sidnal from noise* At low freeuencies? the 
double filter has less phase shift* However? phase shift is 
more severe at freeuencies above 0*5 radians* 

D* PERFORMANCE COMPARISON 



The basic scalar filter and the improved (double) filter 
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were compared usind the data illustrated in fidure 3. A 
series of independent normal random numbers with zero mean 
and variance ten were denerateds- to simulate measurement 
noise* They were added to a nominal function which was a 
combination of steps and a ramp* The ramp function rises 1 
unit each measurement interval r which corresponds to one 
standard deviation in process noise so at least durind the 
ramp function ? the data corresponds to the filter desidn 
variance ratio Q/R ~ 0.1* A spectral analysis of the data is 
presented in fidure 9. The sidnal and noise are presented 
separately and in combination* Note that the noise-only scale 
is expanded* The noise spectrum is irredulari- but overall 
Quite flat* The sidnal consists mostly of very low 
f reQuencies X hut also has some hidh freauencies* This would 
be expected? since step and ramp functions reauire very hidh 
freauencies in their Fourier expansion* The hidh-f reauency 
sidnal is submersed in noise* The data? of course? does not 
fulfil the assumptions from which the Kalman filter is 
derived* However? the real world seldom does either* We are 
lookind for robustness* 

The filters were first tested on the sidnal alone* The 
results are presented in fidure 10* It can be seen that 
neither filter can respond instantaneously to the 
discontinuities in the function? since hidh freauencies are 
attenuated* Both filters lad after discontinuities and 
durind the ramp rise* This is a conseauence of the non- 
symmetric nature of the filters and illustrates the phase 
lad* Note that the double filter performs a bit better? 
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is a pecul isri t'j 



reniainin^ closer to the data throughout <it 
of the plotting routine that only one symbol will be plotted 
when ever data points coincide* The ref o re r data points that 

do not appear should be regarded as occurring simultaneously 
with the ones that do appear)* 

The filters were then tested on the noisy data* The 
results are presented in fisJures 11 and 12* The latter plot 
has the data suppressed so that the scale can be increased 
and more resolution obtained* 

The same trends can be observed as were previously* The 
double filter lads less durinS the trend and transitions* 



The double fi 
closely^ but 
basic filter* 
filter was 4*2 
filters' which 
modification * 
double filter 
remembers less 
even thoudh t 



Iter appears to follow the noise a bit more 
overall it follows the siSnal better than the 
The average variance between the siSnal and the 
0 for the double filter and 5.10 for the basic 
was an 1S% improvement for this simple 
The improvement is due to the fact that the 
weights more recent data more heavily y and 
of the past than the scalar Kalman filter? 
he weisht on the present observation is the 



same * 



This simple 
reader with poss 
any tool? the 
indiscriminant ly 
HamminS C ref *711 



experiment is only intended to acauaint 
ible improvements to the Kalman filter* 
Kalman filter should not be 3pp 
* The interested reader is referred 
for the basic principles of digital fi 
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KIGURE 11 
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E. THE TRANSIENT CASE 



The sophisticated reader has no doubt noticed that the 
steady-state scalar Kalman filter is eoui valent to the 
Bo>c-JenkiriS IMA(Oi»l!'l) model Cref.llll and to Brown's 
exponential smoothing model Cref.103* Zehna criticised the 
exponential smoothing model Cref.133/ notind the bias would 
occur if the steady-state model was applied with an 
inappropr iate prior estimate. Bessler and Zehna Cref.l43 
developed a dain schedule which they call finite exponential 
smoothind. Their formula for dain is 

a ( t ) = a/ ( 1-b^ ) 

where a is the steady state dain)> b = (l-a)y and a(t) is the 
dain schedule as a function of time. It is similar to the 
Kalman dain schedule if the initial Kalman dain K(0) is 
chosen as one. In both modelsi- no prior estimate is 
reouired. The weidht on the first observation is one. A 
comparison of the two models is illustrated in fidure 13? for 
a steady-state dain of 0.2 and an initial dain of 1. The 
Kalman dain was calculated accordind to the recursive formula 
in section I I I. A. 

The Kalman filter dain converses faster? althoudh the 
difference is not dreat. The scalar Kalman filter possesses 
two other advantades over the finite exponential smoothind 
technieue. First? if a dood prior estimate does exist? the 
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Kalman filter allows the initial Jiain to be chosen as less 
than onej and its value will he determined by the assumed 
covariance of the nrior estimate. AIsok the Kalman model 
forces the analyst to at least think about the concepts of 
measurement noise and process noisev and to estimate the 
noise variance ratio Q/R. 

As noted beforey a freauency analysis of the transient 
case is not appropriate. However y it can be thought of as a 
case of transient bandwidth. The filter is initially set to a 
hisfher ^ain than steady-state Sain. If there is no Prior 
estimate availahley the initial Sain is one? and the filter 
is initially an infinite-bandwidth or all-pass filter. As 
data are acauired Sain drops and the bandwidth narrows until 
steady-state conditions are achieved. The concept of 
transient bandwidth is important to the subject of adaptive 
filterinSy to which we will return. 



F. HIGHER ORDER FILTERS 



The main beauty of the Kalman filter is not in its 
statistically unbiased method of calculatins Sainy but in its 
powerful matrix formulationy which allows it to be applied as 
a multi-dimensional model incorporat inS any order of 
differencinS desired. As the state space is increasedy it 
Quickly becomes impossible to analyse the filter 
analytical ly . HiSh-order multi-dimension filters can also 
easily exceed the capacity of present digital computers for 
real-time oppI icat ions . Fortunatelyy it has been found that 
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the state space can he reduced and dimensions decoupled with 
very little degradation in the overall accuracy of the state 
estimate Crefs. -4 and 53. For example!- if a 12-state model 
can be reduced to 9 states and can be adeeuately represented 
bvi three 3-state models» the matrix calculations can be 
considerably simplified and speeded up. 

We will examine a second-order (first difference) 
filters' which can be used to estimate trendi- or velocity. We 
will use the latter term. Position and velocity are to be 
estimated based only on successive measurements of position. 
The state transition and the observation matrices are 



The covariance matrices l.r P> and Q arey of coursey 2 by 
2 matrices. The state vector has two elements y velocity and 
positioriy while the measurement vector has only position. 
The measurement error R is a 1 by 1 matrix which we will 
vary. We have chosen Q as 



arduinS that any process noise will be contained entirely in 
velocity. That isy there can be no random motion that is not 
caused by a random velocity. Randomness of velocity will feed 
into position throudh the state transition matrix. 




H = Cl 0.1 
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Even in thi‘5 



siiiiFle cssey 



solvind snelytical ly for 



istesdy-atate dain in terms of R and Q reouires solvind a 
system of ^th order polynomial eauations. We will opt 
instead for a computer solution* The reader may continue to 
think in terms of the noise variance ratio y where R will take 
on the values ly lOy 100 and Q will remain constant as above* 
Since there is only one non-zero term in the Q matrixy we may 
think of the noise variance ratio as the scalar Quantity 
Q(252)/R* The resulting steady-state dains are 



Noise Variance 
Ratio 



Position 

Gain 



Velocity 

Gain 



1*0 



*769 



*•481 



0*1 *553 *211 

0*01 *362 *080 



As would be expectedy the position dain is much hidher 
than that of a scalar filter at an eaui valent noise variance 
ratio y because the process variation now applies to velocity 
rather than position* The velocity dain is considerably less 
than the position dainy since the velocitvi is not measured 
directly but must be estimated from successive measurements 
of position* The impulse-response function of the medium-dain 
filter (noise variance ratio 0*1) is presented in fidure 14* 
The amplitude response of these three filters is 
compared in fidures 15 and 16* The most strikind feature is 
the ampl if ication which occurs at a specific freouency in the 
position freouency resP'Onse* This imP'lies that the filter is 
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most sensitive to motion in 3 r-^articuler freeuency ran'ie. 
Thusy the nstur<3l freouency of the system to he observed/ if 
it is known/ is s significant design parameter. Adain/ we 
observe the attenuation of hidh-f reeuency noise/ althoudh a 
sidnificant amount still remains at the maximum freouency 
(note that fidure 16 does not include the oridin) . 

The freouency response of velocity shows a reduction in 
amplitude at low freouency. The amplitude at zero freouency 
is eoual to the velocity dain. This is far from ideal 
performance for a differentiator/ which should have an 
amplitude response of zero at zero freouency/ with a slope of 
one UP to the cutoff freouency Cref.7T. The differentiator 
is/ however/ reasonably effective at reducind the amplitude 
of hi dh-f reouency coiTiPonents . 

The phase shift of the filters adain shows increasind 
Phase lad as dain is decreased. The overall effect is 
similar to the scalar filter/ 'and is otherwise unremarkable. 
Therefore/ plots are not included. 

The data of fidure 8 was tested on the lowest-dain 
velocity filter. Note that the spectral content of the data 
(fidure 9> is Quite low-f reouency / and that the bandwidth of 
the lowest dain velocity filter is Quite wide/ and indeed is 
hidher than that of our scalar filter. So it midht be 
expected that the velocity filter would have some trouble 
with the d 

The velocity filter performance on the nominal function 
only is presented in fidure 17. The filter overshoots Quite 
badly at the discontinuous steps/ which/ of course/ are an 
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impulse in velocity. The overshoot is less at the start and 
stop of the trend. It does settle down and track the trend 
without lad> which is an improvement in performance over the 
scalar filter. It should he noted that with hisfher dairu' the 
filter would track the nominal data better/ while with lower 
Sain/ the overshoots would be more severe. 

The performance of the velocity filter on the data is 
illustrated in fidure 18. As expected/ the filter tends to 
follow the noise too much. However/ it does follow the 
discontinuities much more euickly than the scalar filter. 
This points out the fact that the hidher-order filter is more 
effective as a maneuver detector but it is less suitable for 
smoothing very noisy data. This adain illustrates the 
concept of bandwidth/ which is ouite hidh even in the 
low-dain velocity filter. 

In retrospect/ the decision to choose Q<1/1) as zero may 
not have been wise. Allowind some process noise in position/ 
exclusive of velocity/ could well have some smoothind effect 
on the velocity estimate/ which would result in smoother 
one-period ahead p'redictions . This could smooth the 
operation of the filter a bit. The possible combinations of 
filter parameters/ even for this simple filter/ are Quite 
numerous . 

The freauency response of a second-difference 
(acceleration) filter was also determined for comparison. 
The results are presented in fidure 19. The Q matrix was 
zero except for Q(3/3)/ which was one. R was chosen as 10/ 
resultind in a nominal noise variance ratio of 0.1. Adain/ 
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FREaUEMCY RESPONSE OF ACCELERATION FILTER 
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IV 



ESTIMATION y SMOOTHING y AND PREDICTION 



A. ESTIMATION 

So fary we have been concerned only with estiniation of 
the present state. A filter designed to provide such an 
estimate cannot be syinmetricy because it can put no weidht on 
future observations. Thusy phase lad is inevitabley and is 
one of the parameters that should he considered in the desidn 
process . 



B. SMOOTHING 

Smoothind is the use of a filter to provide an estimate 
of p-ast states. Such a filter can be made symmetricy which 
completely eliminates the phase lad. Non-recursive smoothind 
filters cause a loss of N data points at each end of the 
datay where the span of the filter is (-Ny N) . 

The Kalman filter can be used as a smoother by simply 
runnind the forward estimate throudh the filter in the 
opposite direction. The impulse response function of the 
scalar filter was 



d ( t ) = ab^ 



It can he shown (appendi.M A) that the impulse response 
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function of the smoother (forward and backward filters 
combined) is 



G<t) = ab*/< 1+h) 



which is Just the convolution 



3(t) (D 



The Kalman filter is able to provide an estimate 
throughout the span of the data. No data is lost at either 
end. However? due to transient effects? the data near either 
end is subject to phase shift and some increase in Sain. The 
filter is necessarily not symmetric near each end of the data 
span . 

Gelb Lref.'43 includes a complete discussion of 
fixed-point? fixed-las? and fixed-interval smoothinS. We 
will restrict our attention to the scalar? fixed-interval? 
steady-state case? isnorins the end effects. 

The scalar Kalman filter of section III.D (noise 
variance ratio of 0.1? Sain of 0.27) was used as a smoother 
on the data of fiSure 8. The results are presented in 
fiSure 20. As compared to the one-pass performance as 
illustrated in fisure 12? the smoothed data shows phase laS 
removed and peaks in the oscillations reduced. However? the 
smoother has less ability to follow the discontinuities in 
the nominal function. The removal of the phase lad is 
characteristic of any symmetric filter. However? the reduced 
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KAL.MAH SMOOTHliK FtisFOK.' 
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ability to follow the discontinuities is the result of 
reduced effective dsin. Since we have convolved the filter 
weights? we have souared the amplitude of the freouency 
response. The weidht on the data point at (t=0) is reduced. 
The effective sain was 0.27 for the forward filter» and 0.156 
for the smoothinS filter. This reduction in effective Sain 
is not addressed in the literature on the Kalman filtery and 
it is unclear how smoother Sain should he chosen in relation 
to the noise variance ratio Q/R. 

C. A COMPARISON OF TWO SMOOTHERS 

The Kalman smoother of the previous section was compared 
with a Gaussian smoother to illustrate some desiSn options 
and procedures. The Gaussian smoother was chosen from amonS 
a huse variety of data windows because it has Sood smoothinS 
propertiesy and because it is particularly easy to desisn. 
Inte rest indly y preliminary experiments showed repetitive 
applications of a Kalman filter to result in an approximately 
Gaussian filter weisht distribution. A comparison of the 
Gaussian smoother to a variety of other windows is contained 
in Harris Cref.153. 

The Gaussian smoother is a symmetric filter with the 
weights chosen according to a discretized and truncated 
normal distribution. The formula is 

,s(t) = K exp ( -t^/2c*) 
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where t is an inteser on the ranSe (--Ny N) and K is chosen 
such that 



^ si(t) = 1 

The ease of desi.sn comes from the observation that the 
Fourier transform of the continuous Normal distribution is 
also a Normal distribution with scale parameter (variance) 
eaual to As lons:^ as (C3^ > 2) and truncation is not more 
severe than |N| > 2^y a reasonable approximation of the 
freauency response for the Gaussian digital filter is 

G<v) Si gkp(~<t*v^/2) 

The scale parameter was chosen such that the f reeuency 
response was eaual at (v - 0.5)* SkippinsJ the algebraic 
details!» this reauired (6"= 3»11>. The Gaussian smoother was 
truncated to (N = 7) > resulting in a filter span of 15 data 
points. The freauencu response of both filters is presented 
in fiaure 21? and the filter weiahtind coefficients are 
presented in fiaure 22. 

Since we truncated the Gaussian smoother? we would 
expect some ripples in the tail of the freauency response? 
which are Just barely visible in fiaure 21. The Gaussian 
filter has a sharper transition band? and is auite effective 
in blockina hiah freauencies. As compared to the Kalman 
smoother? the Gaussian filter weiahts the present data point 
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lessy nearby data points more and farther data points less. 

The perf oriTiance of the Gaussian smoother on the data of 
fidure 8 is presented in fidure 23. A comparison of the 
Kalman and Gaussian smoothers is presented in fidures 2-4 and 
25. In fidure 24? the smoothers are applied only to the 
nominal function. It can be seen that the Gaussian smoother 
followed the discontinuities and corners of the nominal 
function better than the Kalman smoother. However? when the 
smoothers were applied to the noisy data? the results were 
less clear (fidure 25). The Gaussian smoother adain followed 
the nominal function a bit better? but it also followed 
low-f reauency components of the noise a bit more? tendind to 
emphasize cyclic effects that aren't really there. The mean 
seuare difference between the smoothed estimate and the 
nominal function were very similar? 1.81 for the Gaussian 
smoother and 1.84 for the Kalman smoother. Thus? the Kalman 
filter seems ouite effective when used as a smoother. The 
reader is reminded that the mean seuare difference for the 
scalar Kalman filter was 5.1? which clearly indicates the 
superiority of smoothind over filterind. 



D. PREDICTION 



Prediction is difficult. Recall that a stochastic 
process is an ensemble of possible paths? while data is the 
manif estation of one member of that ensemble. What could 
have happened did? but what can happen isn't necessarily 
doind to. Prediction can be thoudht of as filterind without 
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SCALE FACTOR FOR ORDIMATE 



me<3sureiT»Gn bs y <3nd the estinnste is projected fov'ward in time 
through the state transition model ♦ Filtered or smoothed 
estimates mas be ouite accuratey even if the state transition 
model is not* However? .Sood predictions are heavily 
dependent on an accurate state transition model* 

The his^her-order Kalman filters are polynomial models* 
Hamming has pointed out that polynomial models are poor 
predictors? since the estimate tends to veer off to plus or 
minus infinity as soon as the model is released from the data 
Cref*73* There is no reouirement to use the same model for 
prediction as for filterinS* For example? it mi<^ht make sense 
to track a target with an acceleration filter? hut to compute 
fire control information based on a constant-veloci ty model? 
since tarSet acceleration is Generally assumed random with 
:::ero mean* Similarly? the economist may desire to filter 
data with a hish-order model? but make predictions based on 
constant trend* Clark Cref*53 discusses a somewhat more 
sophist icated method due to SinSer? in which the model decays 
from an acceleration predictor to a constant-veloci ty 
predictor as prediction time increases* Such technioues are 
heuristic in nature? but can prove valuable to the innovative 
analyst * 

An interesting example of the above concept can be found 
in Box and Jenkins Cref*113* They compared a Quadratic 
forecast due to Brown II ref *103 with their own IMA(0?1?1) 
model with a Sain of 0*9* The latter model is eouivalent to 
the steady-state scalar Kalman filter* The data used for the 
comparison was a time series of IBM stock prices* Box and 
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Jenkins observed thaty 



while the eusd ratio uiodel mi^ht well 



be used to fit the data? its perf ormance as a predictor was 
clearly inferior to the simpler IMA(Oylyl) model* This is 
not surprisin.dy since it has lon£:^ been su£^<dested that stocM. 
prices behave as a random walk? and that the best forecast of 
stock pricey at least in the short runy is the present price 
Cref*llZI* Note that the foresoins implies that the dain 
should be set to l*Oy which corresponds to no filtering at 
all* Therefore y Box and Jenkins apparently found that some 
filterini^ of the data was approp^r iate y even thous^h the ^^ain 
they used was Quite hii^h* 
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V. SOME REFINEMENTS^ EXTENSIONS AND ALTERNATIVES 



A. ADAPTIVE FILTERING 



In the linear Kalman filter? the ^ain is completely 
independent of the data* Clearly? this will result in major 
errors if the dain is chosen inapprop r iately or if the data 
statistics chande* If the dain is too low? the filter lads 
badly* In the extreme case? which occurs if the process 
noise covariance Q is much too low? the filter pays much too 
much attention to the past and diverdes from the data* On 
the other hand? if the dain is too hidh? the filter pays too 
much attention to the data and the state estimate contains 
noise* If the filter is a polynomial model and is to be used 
as a predictor? the resultind errors will be spectacular* 

The solution is simple in concept but can be difficult 
to implement* One simply sets the steady-state dain as low 
as appropriate for the stable process beind estimated* In 
tardet trackind? the dain would be set to track an airplane 
flyind a straidht path* A "maneuver detector" or "trend 
detector" is incorporated? which is nothind more than a 
recursive statistical test applied to the residuals to 
determine whether or not they come from a zero-mean 
distribution* If not? the bandwidth is dradually widened 
(dain is increased) until the residuals pass the zero-mean 
test* Then? the dain is allowed to decrease toward the 
stable? steady-state value* 
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Further details and some novel approaches are discussed 



in Clark C reft 5 3 ♦ Two examples taken from Clark are 
illustrated in fisure 26t The conceptual adaptive filter 
discussed above reauires time to detect the maneuver or 
trendy adapt to ity and reconverj^e to a stable ^ain setting ♦ 
Durin^ii this time? the state estimate is less accurate? and 
the time reauired mavj be unacceptably lonS for some 
appl ications ♦ 

Clark proposes a dual-bandwidth adaptive filter to speed 
adaptation* The process is simultaneously tracked by a 
nar row-band and. a wide-band filter* If a maneuver or trend 
is detected? the state estimate of the wide-band filter is 
fed into the narrow-band filter* Ideally? this would allow 
the narrow-band filter to Jump immediately to the current 
(unbiased) estimate of the wide-band filter* In practice? 
Clark found that some widenin<3 of the bandwidth of the 
narrow-band filter was also reauired* 

Voluminous literature exists on the subject? much of it 
very difficult to read* Clark Cref*53 incorporates a 
particularly lucid account of stability problems encountered? 
methods of reducin<^ the cost of false detection of bias? 
analytical methods of determining filter parameters? and 
experimental results* Although Clark's filter was designed 
to track and predict the position of airborne targets? the 
methods discussed are adaptable to the filtering of economic 
time series or virtually any other stochastic process* 
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ACCELEl^ATION ACCELERATIOH 




Conceptual Adaptive Filter 



TIME 




TIUE 



77 



B 



NON-LINEAR FILTERING 



In the non-linesr Kalman filter^ one or more of the 
matrices Q? Ry H? or J are allowed to vary with time* Since 
this results in a tiine-variat ion of the ^ain matrix/ 
GLiantitative analysis in the f reeuency domain is no longer 
appropriate* However? it is well to keep the concepts in mind 
in order to siain additional insis^ht* There are two basic 
types of non-1 inear i t ies that may arise? non-linear 
measurements and non-1 inear dynamics ♦ 

1 ♦ Non-Linear Measurements 

Non-linear measurements arise when observations are 
made in one coordinate system and the model reouires that the 
state he estimated in another coordinate system* In this 
case? the matrices R and H are time-varyins functions of the 

I 

coordinate transf ormation y and do depend on the data? in the 
sense that they depend on the location of the data within the 
coordinate system* This type of non-linearity is often easy 
to handle* 

For the best example of non-linear measurements we must 
return to the tarSet-t racking model* Fire control systems 
generally track in a^rimuth? elevation/ and ranSe* However/ 
the model is a polynomial in Cartesian coordinates/ but not 
in polar coordinates * Airplanes often fly a straight path/ 
but seldom/ if ever/ fly a constant bearinS or ranSe with 
respect to the radar observer*- In this model? the 
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non-1 inearity 



can be reduced by considering the Cartesian 



iiieasurement error as a linear transformat ion of the polar 
measurement error* If the polar measurement error is 
Gaussian i> the Cartesian measurement error is very nearly 
Gaussian with covariance matrix R a function of the 
coordinate transformation* 

The Cartesian R matrix will not he diadonalf even if the 
polar R matrix is* Howeveri> Clark Cref*5J has found that 
settinS the off-diadonal terms of the Cartesian R matrix to 
zero did not appreciably dedrade filter performance* In this 
ways' he was able to decouple a nine-state filter into three 
three-state filters* 

If the measurement non-linearity is too severe » it may 
not be reasonable to assume that the noise is Gaussian* 
However^ limited experiments performed on data with 
non-Gaussian noise (an exponential distribution was used) 
showed that the Kalman smoother and the Gaussian smoother 
were Quite robust as lond as the dain was not hidh* This 
seems to be a conseauence of the Central Limit Theorems- since 
low dain implies a linear combination of a fairly larde 
number of data points* It should he noted that a filter 
desidned to handle this situation is still 1 inear s- althoudh 
the Gaussian assumption is violated* 

2 * Non-Linear Dynamics 

Non-linear dynamics are considerably harder to handle 
than non-linear measurements* This is unfortunates' since the 
areas of potential application are numerous* Non-linear 
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dynaiiiics occur when the Q or ^ matrices depend on the 
previous history of the process. 

As 3 simple example? consider the multiple regression 

model 



y(t) = 3 >J,<t) + h x^(t) + «(t-l) 



where y(t-l) corresponds to the intercept term. In this 
model? we wish to estimate the dependent variable y(t). To do 
so? we need to estimate not only the independent variables 
K,(t) and ;<j(t)? but also the regression coefficients a and b. 
Let us assume we can measure y(D? X|(t)? and x^(t)? but not a 
and b. Assuming a first-order system? the state update 
eauation is 
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or even 



C 1 3(t;t)/2 b(t!t)/2 Kj(t:t)/2 >(^(t,'t)/2 1 

and it obviously changes at every iteration* It is at this 
point that filter desi<=in becomes an art* 

Mote that the independent variables and the reif^ression 
coefficients are assumed here to be first-order Gauss-iiarkov 
processes* Increasingly hidh orders would multiply the state 
space * 

Several experiments were run usiniii a second-order model 
similar to the above on the Box- Jenkins Href* 1 ID series M 
data (sales data with leadinS indicator)* Quantitative 
results are not presented? because the Box-Jenkins data did 
not include sufficient forecast estimates for comparison? 
some "cheating' was done because the Box-Jenkins parameters 

I 

were used in filter design? and it never became clear exactly 
what parameters were appropri<tte for the R and Q matrices* 
However? some Qualitative comments are appropriate* The 
model did work* Some instability was noted in the regression 
parameters* It became obvious that the s^ain on the 
regression parameters must be set very low in comparison to 
the slain on the independent variables? in ordt?r to keep the 
re<^ression parameter estimates from varying faster than the 
estimates of the independent variables* This implies 
choosind small values for the noise variance of the 
rei^ression parameters* Also? bvj keepinS the Sain fairly low 
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on the lesdiniil indicators it was possible to induce a Phase 



la^ that approKijTiately cancelled the lead* 

The reoui rement to keep ^ain low in order to improve 
stability is evidently a conseeuence of the increased decrees 
of freedom* The more parameters to be estimated^ the more 
de<^rees of freedom in the model* His(h <3ain is analaSous to 
relatively few data points bein^ used in res^ression* The 
more variables we introduce into the models the less dain we 
are able to use* 

It is indeed unfortunate that multiple res^ression is a 
non-linear problem when cast in a filter model* It would be 
useful to have a multiple regression model for which more 
recent observations were wei«3hted more heavily than older 
ones in determining the regression parameters* No doubt the 
innovative analyst could develop one to fit the specific 
situation* However s clearly-defined technioues with 
demonstrated results are not available to the pract it ioner * 
The experts all have their favorite methods s and much of the 
literature is difficult to read* There is clearly a need for 
additional research in this area* 

C * NON-PARAMETRIC FILTERING 

We close our discussion with an interestin.^ alternative 
to conventional di^Jital filterin.d technioues* There are 
those who are bothered by the usual distributional 
assumptions made in any application of parametric statistics* 
An extensive literature has developed in the field of 
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non-paraiTietric statisticsj^ which is based on the principle 



that distributional assumptions are avoidedy or at least 
weakened^ A stronsi point of non-paramet ric statistics is 
relative insensitivity to extreme outliers^ Howevery little 
pro.^!ress has been made in the non-paramet ric analysis of time 

series* An exception may be found in Tukey Cref*163y which 
is presented in a highly intuitive mannery with little or no 
theoret ical background ♦ 

One simple idea advanced by Tukey is that of median 
smoothirii^* The smoothed estimate is based on the median of 
several adjacent data pointsy rather than on a weighted 
linear combination* The result is obviously a series of 
steps y since adjacent data points will often have the same 
median* Tukey su^dests several methods to restore some 
curvature in the estimate* These will not be developed here* 
Tukey 's methods would be relatively hard to mechanize on a 
computer y because the methodolosiy reouires extensive logical 
rules * 

Tukey 's methods could he extended to real-time filtering 
problems by developing a non-parametr ic analos to the 
recursive dii^ital filter* Recall that the recursive digital 
filter consists of a weighted linear combination of recent 
data points added to a weighted linear combination of recent 
estimates* The non-parametr ic filter estimate would simply 
be the median of several recent data points and several 
recent states* The idea is intuitively appealing and should 
be the subject of future research* Discussion here will be 
limited to some of the more obvious traps that await the 
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unw3 ry . 



If such 3 filter were to he desidned» the 
impulse- response function would he inesnindless > hecsuse the 
iTiedien estimate would always be zero if the span of the 
filter were Greater than two t Nevertheless^ a median filter 
does have a freauency responser which in fact is a 
particularly nasty one. 

Consider a seven-point median smootherr where the state 
estimate at time t is the median of the measurements made at 
time (t-3) to time (t+3). This is analadous to the 
rectangular (parametric) window discussed in Hammind Cref.73. 
The rectandular window weights all data points within the 
window eaually. The median window obviously does the same. 
As a result? we would expect the freauency response of the 
median window to have severe ripples as does that of the 
rectangular window. We can see intuitively that this is true. 
Since the span of our example median window is 7r the 
freauency response of any freauencvj that is a non-pero 
integer multiple of 1/7 is obviously zero. The freauency 
response at zero freauency is one? since the zero freauency 
implies a constant. The amplitude of the freauency response 
falls off to the first zero? then rises again. Successive 
maxima decrease with increasing freauency? but the freauency 
response is always non-zero except at freauencies that are 
non-zero integer multiples of the reciprocal of the span of 
the window. Thus? the non-paramet ric filter will need to 
incorp-orate some (unspecified) device to improve the 
freauency response. 
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Two other difficulties ere worthy of mention 



F i r s t f 



the seiiiplind distribution of the medien may have a larder 
variance than the samplind distribution of the mean. This 
means that the parametric filter may provide a better 
estimate than the non-parametric filter if the assumptions on 
which the parametric filter is based are at all reasonable. 
Secondf for the non-parametric filter to be useful » the 
median must be a statistic of interest. If it is assumed 
that the distribution is symmetrici' the median and mean aref 
of courser eoual. If the samplind distribution is skewedr the 
mean cannot be deduced from the median unless strict 
parametric assumptions are imposedr which of courser override 
the Justification for the non-parametric filter in the first 
place. The idea is nevertheless intriduindr and should be 



explored further. 
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SOME APPLICATIONS 



The Ksluian filter has been applied to Operations Research 
and economic problems with var^ins^ decrees of success^ 
hcUhorter C ref *1711 conducted an empirical study of the Kalman 
filter in which he compared it to several other methods of 
time series forecasting* The results were mixed? with no 
method dominating* The Kalman filter compared more favorably 
over a short term forecasting horizon than over a lond term 

one* Its performance was? not surpr isini^ly ? found to be 
degraded if the structural model was seriously mis“**specif ied * 
McWhorter pointed out some of the difficulties encountered in 
building the model* In an economic context? it is often very 
difficult to specify the noise covariance matrices R and Q? 
and even to identify the structure of the state transition 
matrix ^* The assumptions made are often sweeping and 
arbitrary? in contrast to trackini^ applications where the 
noise processes and especially the state transition model are 
relatively well understood* 

A* INVENTORY MANAGEMENT 

The Kalman filter is directly applicable to inventory 
manasement? and if properly desi<^ned? should be superior to 
the finite exponential smoothing model of Bossier and Zehna 
Cref*14H* DowninS? Pike? and Morrison C ref *183 designed a 
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Kalman filter for the inventory control of nuclear material 



The paper is readablev and the filter is well-documented and 
easy to understand* They use the concept of a control 
vector » which has not been mentioned here* An interesting 
peculiarity of the model is that one of the measurements is 
only available once every twenty iterations* The state 
transition matrix is a simple material balance relation which 
is obviously Quite accurate* Such a model could be expected 
to perform Quite well* 

B* ESTIMATING A MEAN FUNCTION 

Although the Kalman filter was derived from an assumption 
of stationarity i» we have seen that it can be Quite powerful 
in separating a time varyind sidnal from noise* The examples 
of section III were all essentially estimates of the 
time-varyind mean function of a stochastic process* The 
example process was Gaussian with a constant variance* The 
variance was the measurement noise and so directly 
influenced the dain* If variance were not constants the 
performance of a non-adaptive filter would be dedraded* If 
the chande in variance was dreat enoudhr an adaptive filter 
would be reauired* 

A dood method of estimatind a time-varyind mean function 
could be applied in numerous areas > such as any sort of 
traffic or flow control problem^ perhaps in Quality control 
of larde-batch or flow manuf actur ind processes? and any 
application where it is desirable to detect a chande in the 
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process* The sensitivity of the filter is directly 
adjustable by the modeller through the noise covariance 
matrices Q and R* 

A particularly useful application would he to the 
estimation of the rate Parameter of a non-t ime-homo^^eneous 
Poisson process* If this can he done accurately? the process 
can be transformed to a stationary one C ref ♦63? which if^reatly 
ejipands the number of analytic tools that can he used* 

The Poisson process is a counting process in continuous 
time? and to attempt to filter a strinsi of interarrival time 
data would violate the samplinsJ theorem* The times of 
arrival are the measurement times? and they are most 
certainly not made at eoually spaced intervals* Instead? the 
filter may be desis^ned to sample a countin*^ process* At 
discrete intervals the filter would count the number of 
arrivals since some arbitrarvj time ori?^in* If the F'rocess 
were to continue for a lonsi time? the time origin mis^ht 
occasionally have to be reset to prevent computer overflow* 
It is easy to see how this samplin.d process could be 
implemented even if the input data were actually arrival 
instants in continuous time* The sampling interval should be 
small enousih that there is low probability that more than one 
arrival would occur durind a diven measurement interval* 
Since the number of arrivals is monotone non-decreasind in 
time? a velocity or trend model would be appropriate* The 
input data would consist of inteders* The state estimates 
would not* The non-inteder estimate of number of arrivals up 
to the current time would not he useful to us* However? the 
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second element of the state vector? the velocity or trend? 
would in fact be the filtered arrival rate estimate. Since 
the process is noisy and non-Gaussian ? a very low 
steady-state Sain is appropriate. 

The time-varyind Poisson process cannot have constant 
variance? since the mean and the variance are eoual. A low 
arrival rate implies hidh variance in the Poisson process? 
which is eauivalent to hidh measurement noise? which reouires 
low dain. A constant-dain filter would therefore be 
relatively more sensitive at low arrival rates than at hidh 
arrival rates. An adaptive filter could be easily desidned 
to use the inverse of the rate estimate as the measurement 
noise variance estimate. Stability midht reouire that the 
adjustment of the measurement noise variance be itself a 
filterind process? in which the incomind variance estimate is 
redarded as data. 

C. MULTIPLE REGRESSION 

If the redression constants are assumed known <or 
computed by other means) the desidn of an appropriate filter 
is Quite st raidhtf orward ? and Quality of estimation is 
related directly to the Quality of the model. Note that the 
velocity filter is simply the redression of velocity on 
position? where the slope parameter is known to be one. If 
the redression coefficients are assumed to vary in time? the 
problem becomes non-linear and is Quite complex. Because of 
the immense applicability of this model? additional 
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developirientsl work is 



indeed a fertile field for future 



research ♦ 

D, SOME DESIGN CONSIDERATIONS 

In applications where the noise covariance matrices R and 
Qj- and the system dynamic model (state transition matrix 
are known or easily estimated^ desidn is straightforward 
and has been successfully accomplished while remainind in the 
time domain* However » in applications where sweepind 
assumptions are reauired» a f reouency -domain analysis could 
be very helpful* Some duidelines are as follows J 

1 * Spectral Analysis of the Data 

A spectral analysis of sample data will show what the 
freauency response of the filter should he* The Fast Fourier 
Transform (FFT) prodram available in most computer libraries 
is denerally easy to use* However? the FFT prodrams 
denerally reauire an exact power of 2 for the number of data 
points* Hammind Cref*?!] points out some pitfalls* Since 
stationarity is assumed? the data should be considered as a 
rotatind cylinder? and if the startind and endind values are 
not similar? a discontinuity will exist in the spectrum* The 
data can be tapered and padded with zeros? but exactly the 
best method to accomplish this is unknown* Several methods 
midht be tried* 

The main virtue of the FFT is its speed* It works well 
on a lond run of data* If the number of data points is small 
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(around a hundred) it midht be effective to find (or write) a 



less efficient.* conventional discrete Fourier transforiri 
prodraiTi* which would not reeuire padded* truncated* or 
tapered data if the starting and endinS values are similar* 
If the FFT prodram used does not reeuire an exact power of 
two for the number of input data points* it would be well to 
find out why not. The prodram may be doind the paddind and 
tsperind itself* and the analyst should be curious as to how. 
The analyst should remember that the spectrum is computed 
from the data* and it is therefore an estimate. If the run 
of data is short* there will be considerable variance in the 
estimate . 

2 . Frequency Analysis of Proposed Models 

The analyst may test the effect of assumptions made 
in desidnind the filter by simply obtainind an impulse 
response of the filter and runnind it throudh an FFT. 
Truncation and taperind is no problem* because the impulse 
response will approach zero with time. The proper impulse 

n 

function is simply a 1 followed by 2 -1 zeros for a filter* 
or a 1 in the middle of 2^ -1 zeros for a smoother* If the 
output of the FFT consists of real and imadinary components* 
it will be necessary to compute amplitude and phase. 

3 . AdJustind the Model 

If the model dynamics seem adequate but the bandwidth 
is wrond* the analyst should by now have some insidht into 
what adjustments to make to the noise covariance matrices to 
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try to iiiiprove thinsis. In 3 model of 3nvj complexity 3t sllf 
there are numerous possible combinations. However > even some 
improvement over the initial assumptions will be beneficial. 
Ule are not lookind for theoretical eledance? we are lookind 
for performance. 

Perhaps the model dynamics obviously call for a trend 
filter or even a chande -of-t rend ( acceleration ) filterf but 
the data is auite noisy. Consideration should be diven to 
lowerind the order of the filter. A very low-dain velocity 
filter will not follow chandes in trend well. A hidhei — dain 
scalar filter may do so more effectively^ althoudh it will 
lad a steady trend. There are many tradeoffs » and we cannot 
achieve perfection. 



4. Testind the Model 
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APPENDIX A. DERIVATIONS 



1. SCALAR KALMAN FILTER 

3 » Recursive Formul3 for Ksliiisn Gain 

The cov3ri3nce extrapolation eoustion 
P( t) = f £( t-1 ) + Q 

reduces in the scalar case to 

P(t) = 2(t-l) + Q 

Since 

K(t) = 2(t) H^r'* 

we may write? for the scalar case? 

2(t-l) = K(t-l) R 

Similarly? since 

T r 

K(t) = PH C HPH + RO 

by reducing to the scalar case and substituting? we may write 

K(t) = K(t-l) R + Q = K(t-l) t Q/R 

K(t-l) R + 0 f R K(t-l) f Q/R + 1 

b . Steady-State Kalman Gain 

RearranSinS the recursive Sain eauation and lettins 
K(t) = K(t-l) = K 

we see that 

K^+ (Q/R)K - Q/R = 0 

By the Quadratic formula? 

K •■= j;:Q + / Q^+ Q_ 

2R / ^ R 



Ue are oblised to take the larSer root? since the smaller 
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root would force the to be negative* We also observe the 



inverse relationship 



z 



Q^/4R^+ Q/R 



(K + Q/2R) 



k'+ (Q/R)K + Q^/4R^ = 0*/4R‘*‘+ Q/R 



K 



..Z 



(Q/R) ( 1-K) 



Q/R = K*/(l-K) 



c ♦ Transient Kalman Gain 

Recall that the Kalman filter reauires a prior state 
estimate X(0) and a prior estimate of covariance P(0). This 
reauirement can be avoided by usinsi K(l) = 1? which allows 
the initial state estimate to be eoual to the first 
measurement* Recall that 



T 

K(t) = I(t) HR 



Since K(l) = 1? then 2(1) - R 



d* Amplitude and Phase of Freouency Response 



The freouency response is 



H( V) 




since 



)b exp(-iv)| < 1 



then 



H(v) = a / IZl-b exp(-iv)J 



The amplitude souared is 



Z 

A ■= H(v) H(-v) = a 



a^/ Cl-b exp(-iv)l Cl-b exp(iv)l 
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I 












By Euler ''s relation 

~ B / Cl f -2b cos V.1 

which msy be written 

A^ - / C ( l~b )* f 2b ( 1 -cos v ) J 

recal 1 in<3 that 

Q/R •■= K^/(l-K) = 3*/b 

we may write the amplitude as 



A 



Q/R 



Q/R + 2 (1-cos 



V ) 



The phase andle is 

Q(v) = arctan CIin(v) / Re(v)] 

where Ini(v) and Re(v) are the imaginary and real parts of 

H(v)» which mavj be written 

H ( V ) = a Cl-b eMp(iv)ZI 

Cl-b e:;p ( -i V ) 2 C 1-b exp(iv)] 

H ( V ) = a ( 1 fb cos y -ib sin v) 
i r~b* “2b cos V 

which allows us to write 

0(v) = arctan C(-b sin v) / (1-b cos v)] 

The andle for maximum p-hase shift occurs when 

d8 ( V ) = - b cos V ~ 0 

dv ~i -f b'*’ -2b cos V 

so that the maximum phase shift 0(v)max occurs when 

V - arccos b 

and has a value of 

0(v)max = arctan -b sin (arccos b) 

1-b cos ( a rccos b ) 

= arctan ( -b / yi-b^ ) 
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IMPULSE RESPONSE FUNCTIONS 



The impulse-resFonse of the sceler Kalman filter is 
S ( t ) = 3b* t = 0 y 1 » 2 ? . » . 



3» Impulse- response of Double Filter 



:^(t) - 3b 0 ab 
3(t) = £a.b‘‘ab‘-'‘ = a" 

a(t) = (t+i) 3*b^ 



IixiPulse-resPonsG 


' of Scalar Kalnian 


Smoother 




t ^ -t 

S ( t ) = ab 0 ab 




aO 


«o 


ao 


d ( t ) = ^ a 


*» tf A. 2 ^ k}*^**" 

b ab =3 ^ b = 


a' b*- £ b'"" 






jb-o 



^(t) = a^b*^/ (1-b^ = a b*/(l+b) 
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